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ABSTRACT 

We perform direct numerical simulations of an externally driven two- 
dimensional magnetohydrodynamic system over extended periods of time 
to simulate the dynamics of a transverse section of a solar coronal loop. 
A stationary and large-scale magnetic forcing was imposed, to model the 
photospheric motions at the magnetic loop footpoints. A turbulent stationary 
regime is reached, which corresponds to energy dissipation rates consistent with 
the heating requirements of coronal loops. 

The temporal behavior of quantities such as the energy dissipation rate 
show clear indications of intermittency, which are exclusively due to the 
strong nonlinearity of the system. We tentatively associate these impulsive 
events of magnetic energy dissipation (from 5 x 10 24 erg to 10 26 erg) to the 
so-called nanoflares. A statistical analysis of these events yields a power law 
distribution as a function of their energies with a negative slope of 1.5, which 
is consistent with those obtained for flare energy distributions reported from 
X-ray observations. 

Subject headings: Sun: flares — MHD — turbulence 
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1. Introduction 

Coronal loops in active regions are likely to be heated by Joule dissipation of highly 
structured electric currents. The formation of small scales in the spatial distribution of 
electric currents is necessary to enhance magnetic energy dissipation and therefore provide 
sufficient heating to the plasma confined in these loops. Various scenarios of how these fine 
scale structures might be generated have been proposed, such as the spontaneous formation 
of tangential discontinuities flParker TWT2, [Parker 1983[ ), the development of an energy 



cascade driven by random footpoint motions on a force-free configuration ( |van Ballegooijen] 
|1986|) , or the direct energy cascade associated to a fully turbulent magnetohydrodynamic 
(MHD) regime (|Heyvaerts fc Priest 1992^ , |G6mez fc Ferro Fontan 1992|) . These rather 



different models share in common the dominant role of nonlinearities in generating fine 
spatial structure. 

In this paper we assume that the dynamics of a coronal loop driven by footpoint 
motions is described by the MHD equations. Since the kinetic (R) and magnetic (S) 
Reynolds numbers in coronal active regions are extremely large (R ~ S ~ 10 10-12 ), we 
expect footpoint motions to drive the loop into a strongly turbulent MHD regime. 

Footpoint motions whose lengthscales are much smaller than the loop length cause 
the coronal plasma to move in planes perpendicular to the axial magnetic field, generating 
a small transverse magnetic field component. In §2 we model this coupling to simulate 
the driving action of footpoint motions on a generic transverse section of the loop. The 
numerical technique used for the integration of the two-dimensional MHD equations is 
described in §3. In §4 we report the energy dissipation rate that we obtain, and a statistical 
analysis of dissipation events is presented in §5. Finally, the relevant results of this paper 
are summarized in 56. 
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2. Forced two-dimensional magnetohydrodynamics 

The dynamics of a coronal loop with a uniform magnetic field B = B Q z, length L and 
transverse section (2nl) x (2ttI), can be modeled by the RMHD equations (|Strauss 1976|) : 

d t a = v A d z ip + [ip,a}+ i]V 2 a (1) 

d t w = v A d z j + [tp, w] - [a,j] + uV 2 w (2) 

where v A = B /y/AiTp is the Alfven speed, v is the kinematic viscosity, r\ is the plasma 
resistivity, ip is the stream function, a is the vector potential, w = — V 2 ^ is the fluid 
vorticity, j = — V 2 a is the electric current density and [u, v ] = z ■ Vu x Vw. For given 
photospheric motions applied at the footpoints (plates z = and z = L) horizontal velocity 
and magnetic field components develop in the interior of the loop, given by v = V x (zip) 
and b = V x (za). 

The RMHD equations can be regarded as describing a set of two-dimensional MHD 
systems stacked along the loop axis and interacting among themselves through the v A d z 
terms. For simplicity, hereafter we study the evolution of a generic two-dimensional slice of 
a loop to which end we model the v A d z terms as external forces (see [Einaudi et al. 1996 
for a similar approach). We assume the vector potential to be independent of z and the 
stream function to interpolate linearly between ip(z = 0) = and ip(z = L) = where 
ty(x,y,t) is the stream function for the photospheric velocity field. These assumptions yield 
v A d z ip = v A ^/L (in Eqn ([[])) and v A d z j = (in Eqn (0)) and correspond to an idealized 
scenario where the magnetic stress distributes uniformly throughout the loop. The 2D 
equations for a generic transverse slice of a loop become, 

dta = [V>, a] + r]V 2 a + f (3) 

d t w = [-0, w) — [a,j] + vV 2 w (4) 

where / = (v A / L)^> . 
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3. Numerical procedures 

We performed numerical simulations of Eqs (H)-(fD on a square box of sides 2tt, with 
periodic boundary conditions. The magnetic vector potential and the stream function 
are expanded in Fourier series. To be able to perform long-time integrations, we worked 
with a moderate resolution version of the code (96 x 96 grid points). The code is of the 
pseudo-spectral type, with 2/3 de-aliasing flCanuto et al. 19881) . The temporal integration 



scheme is a fifth order predictor-corrector, to achieve almost exact energy balance over our 
extended time simulations. 

We turn Eqs into a dimensionless version, choosing I and L as the units 

for transverse and longitudinal distances, and Vq = \/fo as the unit for velocities 
(io = v A u p (l/L), u p : typical photospheric velocity), since the field intensities are determined 
by the forcing strength. The dimensionless dissipation coefficients are vq = vjil vTo) and 
Vo = V7o)- The forcing is constant in time and non zero only in a narrow band in 
fc-space corresponding to 3 < k I < 4. In spite of the narrow forcing and even though 
velocity and magnetic fields are initially zero, non-linear terms quickly populate all the 
modes across the spectrum and a stationary turbulent state is reached. 



4. Energy dissipation rate 

To restore the dimensions to our numerical results, we used typical values for the solar 
corona: L ~ 5 x 10 9 cm, I ~ 10 8 cm, v p ~ 10 5 cm/s, B ~ 100 G, n ~ 5 x 10 9 cm -3 and 
t/ = J)o = 3 x 10~ 2 . 

Fig. 1 shows magnetic and kinetic energy vs time. The behavior of both energies is 
highly intermittent despite the fact that the forcing is constant and coherent. This kind 
of behavior is usually called internal intermittency, to emphasize the fact that the rapid 
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fluctuations are not induced by an external random forcing. Also, note that magnetic 
energy is about one order of magnitude larger than kinetic energy. 

EDITOR: PLACE FIGURE HERE. 

The energy dissipation rate is also a strongly intermittent quantity as shown in Fig. 2. For 
turbulent systems at large Reynolds numbers, the dissipation rate in the stationary regime 
is expected to be independent of the Reynolds number R ( Kolmogorov 194"! ). For the rather 



moderate Reynolds number simulations reported here, a weak (monotonically increasing) 
dependence of the dissipation rate with R still remains (see also Einaudi et al. 19961) . 

EDITOR: PLACE FIGURE g HERE. 

The total dissipation rate for the typical numbers listed above is e ~ 1.6 x 10 24 erg.s -1 . 
We can transform the heating rate into an energy influx from the photosphere, by 
simply dividing by twice the transverse area (because we have two boundaries), i.e. 
T = e/(2(27rZ) 2 ). In Eqn (||) we show the quantitative value of the energy flux as well as 
the explicit dependence with the relevant parameters of the problem, 

t =2x i 6^^ ( n \ k ( B ° y ( u P h ^ (tqHj 1 (5) 

cm 2 s V5 x 10 9 cm" 3 / V100 GJ U0 5 cms-V f k a \\ y! 

V5xl0 9 cm' 

This energy flux compares quite favorably with the heating requirements for active regions, 
which span the range T = 3 x 10 5 — 10 7 erg cm -2 s -1 (|Withbroe fc JNoyes 1977|). 



5. Distribution of nanoflares 



In this section, we associate the peaks of energy dissipation displayed in Fig. 2 to the 
so-called nanoflares ( [Parker 1988|) . We estimate the occurrence rate for these nano-events, 
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i.e. the number of events per unit energy and time P(E) = dN/dE so that, 

rEmax 

K= / dE P(E) (6) 



is the total number of events per unit time and 

(■Emax 

e= dE E P(E) (7) 

/ p 

is the total heating rate (in erg s _1 ) contributed by all events in the energy range 
[E min ; E max ]. A simple inspection to the e(t) time series shown in Fig. 2 indicates 
that these events are in a highly concentrated or piled-up regime, i.e. their event rate 1Z 
multiplied by their typical duration is much larger than unity. At any given time, many 
events are going on simultaneously. 

This piled-up scenario is a serious drawback against performing any kind of statistical 
analysis. To overcome this difficulty we define an event in the following fashion: (1) we set 
a threshold heating rate €q on the time series displayed in Fig. 2, of the order of its time 
average, (2) events are excesses of dissipation which start when e(t) surpasses eo and finish 
when e(t) returns below e . Once a particular threshold is set, we perform a statistical 
analysis of the events, keeping track of their peak values, durations and total energy content. 
The implicit assumption behind our working definition, is that the small fraction of events 
that emerge over the threshold are statistically representative of the whole set. We do not 
make any attempt to prove this assertion, which therefore remains as a working hypothesis. 

Among the interesting results of this statistical analysis, we find a significant correlation 
between the energy and duration of events like E ~ r 2 (|Dmitruk fc Gomez 1997 | ) , which is 



consistent with the correlation reported by |Lee et al. 1993| from hard X-ray observations 
and by [Lu et al. 1993| for an avalanche model of flare occurrences. Perhaps the most 



important result is that the occurrence rate as a function of energy displays a power law 
behavior 

P(E) = A E~ l - 5±0 - 2 (8 
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in the energy range spanned from E min ~ 5 x 10 24 erg to E max ~ 10 26 erg. 

We define the constant A in Eqn ([8]) so that the heating rate computed from Eqn 
(0) matches the total heating rate from our simulation (see Eqn ([5])). Figure 3 shows our 
occurrence rate, displaying the power law behavior indicated in Eqn (|8]). For comparison, 
we also plotted the occurrence rate for transient brightenings derived by |Shimizu 1995 
(slope between 1.5-1.6) from Yohkoh soft X-ray observations and by |Crosby et al. 1993 
(slope 1.53) from SMM hard X-ray data. 



EDITOR: PLACE FIGURE | HERE. 



Also, we obtained the distribution of events as a function of peak fluxes which is a power 
law with slope 1.7 ±0.3. This slope is consistent with the one derived by |Crosby et al. 1993 



(1.68) from X-ray events and somewhat flatter than those reported by [Hudson 199 1| (1.8). 
Note that the slopes derived from X-ray observations assume that the luminosity in X-rays 
is proportional to the dissipated energy. [Porter et al. 1995] criticize this assumption after 



comparing the UV and X-ray emission for several microflares and find that slopes derived 
from X-rays become slightly steeper when corrected for this effect. 

The remarkable correspondence between the rates plotted in Fig. 3, is indicative of the 
presence of a common physical process behind the dissipation events ranging from 10 24 erg 
or less, to 10 33 erg for the largest flares. Since in all cases the index of the power law remains 
smaller than two, Eqn (|7]) implies that the contribution to energy dissipation in a given 
energy range ([E min ; E max \) is dominated by the most energetic events (i.e. E ~ E max ). 
According to this result, the relatively infrequent large energy events contribute more to 
the heating rate than the much more numerous small energy events. This assertion might 
only change if indications of a turn-up of the slope (to an index larger than two) are found 
at low energies ([Hudson 1991| ). 
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6. Discussion and conclusions 

In the present paper we simulate the dynamics of a transverse section of a solar coronal 
loop through an externally driven two-dimensional MHD code. The relevant results of 
this study are summarized as follows: (1) for an external forcing which is narrowband 
in wavenumber, we find that the system becomes strongly turbulent and after about 
twenty photospheric turnover times a stationary turbulent regime is reached, (2) the 
energy dissipation rate obtained for typical footpoint velocities is consistent with the power 
necessary to heat active region loops {T ~ 2 x 10 6 erg cm -2 s _1 ), (3) the energy dissipation 
rate displays a highly intermittent behavior, which is a ubiquitous characteristic of turbulent 
systems, (4) we associate the elementary events that compose this intermittent heating rate 
to the so-called nanoflare events described by |Parker 1988| . A statistical analysis of the 
events performed on a long term numerical simulation (about 200 turnover times) shows a 
power law event rate, going like dn/dE ~ E' 1 ' 5 , which is remarkably consistent with the 
statistics of flare occurrence derived from observations QHudson 1991| , |Crosby et al. 1993| , 



Lee et al. 1993 , Shimizu 1995 ). 
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Fig. 1. — Magnetic and kinetic energy as functions of time for our simulation. 
Fig. 2. — Magnetic energy dissipation rate as a function of time. 

Fig. 3. — Number of events per unit energy and time. In the upper left corner we plotted 
the histogram for our simulation and the thick line shows the best power law fit (slope 
1.5 ± 0.2). Dashed trace correspond to the occurrence rates derived by Shimizu 1995, and 
the dot-dashed trace to the distributions computed by Crosby et al. 1993 
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